#ifdef USE_EIGEN3
#include <Element.d/NonLinearity.d/NLPentahedral.h>

template<>
const double SolidElementTemplate<Wedge6LagrangePolynomialShapeFunction,6,6>
::nodeRefCoords[6][3] = 
  {{0.0,0.0,-1.0},{1.0,0.0,-1.0},{0.0,1.0,-1.0},
   {0.0,0.0, 1.0},{1.0,0.0, 1.0},{0.0,1.0, 1.0}};

template<>
const double SolidElementTemplate<Wedge15LagrangePolynomialShapeFunction,15,9>
::nodeRefCoords[15][3] = 
  {{  0., 0., -1.},{  1.,  0., -1.},{ 0.,  1., -1.},
   {  0., 0.,  1.},{  1.,  0.,  1.},{ 0.,  1.,  1.},
   { 0.5, 0., -1.},{ 0.5, 0.5, -1.},{ 0., 0.5, -1.},
   { 0.5, 0.,  1.},{ 0.5, 0.5,  1.},{ 0., 0.5,  1.},
   {  0., 0.,  0.},{  1.,  0.,  0.},{ 0.,  1.,  0.}};

template<>
const double SolidElementTemplate<Wedge26LagrangePolynomialShapeFunction,26,18>
::nodeRefCoords[26][3] = 
  {{  0. ,  0. ,  -1.},{  1. ,  0. ,  -1.},{  0. ,  1. ,  -1.},
   {  0. ,  0. ,   1.},{  1. ,  0. ,   1.},{  0. ,  1. ,   1.},
   {1./3.,  0. ,  -1.},{2./3.,  0. ,  -1.},{2./3.,1./3.,  -1.},
   {1./3.,2./3.,  -1.},{  0. ,2./3.,  -1.},{  0. ,1./3.,  -1.},
   {1./3.,  0. ,   1.},{2./3.,  0. ,   1.},{2./3.,1./3.,   1.},
   {1./3.,2./3.,   1.},{  0. ,2./3.,   1.},{  0. ,1./3.,   1.},
   {  0. ,  0. ,-1./3},{  0. ,  0. , 1./3},{  1. ,  0. ,-1./3},
   {  1. ,  0. , 1./3},{  0. ,  1. ,-1./3},{  0. ,  1. , 1./3},
   {1./3., 1./3. ,-1.},{1./3., 1./3. , 1.}};

template<>
void
SolidElementTemplate<Wedge6LagrangePolynomialShapeFunction,6,6>
::getGaussPointAndWeight(int n, double *point, double &weight)
{
  double weight3d6[6] = { 0.1666666666666667, 0.1666666666666667, 0.1666666666666667,
                          0.1666666666666667, 0.1666666666666667, 0.1666666666666667 };
  double gauss3d6[6][3] = { {0.666666666666667, 0.166666666666667, -0.577350269189626},
                            {0.166666666666667, 0.666666666666667, -0.577350269189626},
                            {0.166666666666667, 0.166666666666667, -0.577350269189626},
                            {0.666666666666667, 0.166666666666667, 0.577350269189626},
                            {0.166666666666667, 0.666666666666667, 0.577350269189626},
                            {0.166666666666667, 0.166666666666667, 0.577350269189626} };

  point[0] = gauss3d6[n][0];
  point[1] = gauss3d6[n][1];
  point[2] = gauss3d6[n][2];
  weight = weight3d6[n];
}

template<>
void
SolidElementTemplate<Wedge15LagrangePolynomialShapeFunction,15,9>
::getGaussPointAndWeight(int n, double *point, double &weight)
{
  double weight3d8[9] = { 0.092592592592593,0.092592592592593,0.092592592592593,
                          0.148148148148148,0.148148148148148,0.148148148148148,
                          0.092592592592593,0.092592592592593,0.092592592592593};
  double gauss3d8[9][3] = { {0.166666666666667,0.166666666666667,-0.774596669241483},
                            {0.666666666666667,0.166666666666667,-0.774596669241483},
                            {0.166666666666667,0.666666666666667,-0.774596669241483},
                            {0.166666666666667,0.166666666666667,0.},
                            {0.666666666666667,0.166666666666667,0.},
                            {0.166666666666667,0.666666666666667,0.},
                            {0.166666666666667,0.166666666666667,0.774596669241483},
                            {0.666666666666667,0.166666666666667,0.774596669241483},
                            {0.166666666666667,0.666666666666667,0.774596669241483} };


  point[0] = gauss3d8[n][0];
  point[1] = gauss3d8[n][1];
  point[2] = gauss3d8[n][2];
  weight = weight3d8[n];
}

template<>
void
SolidElementTemplate<Wedge26LagrangePolynomialShapeFunction,26,18>
::getGaussPointAndWeight(int n, double *point, double &weight)
{
  double gauss3d9[18][3] = {
     { 0.166666666666667, 0.166666666666667, -0.774596669241483 },
     { 0.166666666666667, 0.666666666666667, -0.774596669241483 },
     { 0.666666666666667, 0.166666666666667, -0.774596669241483 },
     { 0.000000000000000, 0.500000000000000, -0.774596669241483 },
     { 0.500000000000000, 0.000000000000000, -0.774596669241483 },
     { 0.500000000000000, 0.500000000000000, -0.774596669241483 },
     { 0.166666666666667, 0.166666666666667, 0.0 },
     { 0.166666666666667, 0.666666666666667, 0.0 },
     { 0.666666666666667, 0.166666666666667, 0.0 },
     { 0.000000000000000, 0.500000000000000, 0.0 },
     { 0.500000000000000, 0.000000000000000, 0.0 },
     { 0.500000000000000, 0.500000000000000, 0.0 },
     { 0.166666666666667, 0.166666666666667, 0.774596669241483 },
     { 0.166666666666667, 0.666666666666667, 0.774596669241483 },
     { 0.666666666666667, 0.166666666666667, 0.774596669241483 },
     { 0.000000000000000, 0.500000000000000, 0.774596669241483 },
     { 0.500000000000000, 0.000000000000000, 0.774596669241483 },
     { 0.500000000000000, 0.500000000000000, 0.774596669241483 }
  };

  double weight3d9[18] = {
     0.083333333333333, 0.083333333333333, 0.083333333333333,
     0.009259259259259, 0.009259259259259, 0.009259259259259,
     0.133333333333333, 0.133333333333333, 0.133333333333333,
     0.014814814814815, 0.014814814814815, 0.014814814814815,
     0.083333333333333, 0.083333333333333, 0.083333333333333,
     0.009259259259259, 0.009259259259259, 0.009259259259259
  };

  point[0] = gauss3d9[n][0];
  point[1] = gauss3d9[n][1];
  point[2] = gauss3d9[n][2];
  weight = weight3d9[n];
}
#endif
